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Abstract 

The elliptic Gaudin model describes completely anisotropic spin systems with long range inter¬ 
actions. The model was proven to be quantum integrable by Gaudin and latter the exact solution 
was found by means of the algebraic Bethe ansatz. In spite of the appealing properties of the 
model, it has not yet been applied to any physical problem. We here generalize the exact solution 
to systems with arbitrary spins, and study numerically the behavior of the Bethe roots for a system 
with three different spins. Then, we propose an integrable anisotropic central spin model that we 
study numerically for very large systems. 
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I. INTRODUCTION 


In 1976 Michel Gaudin raoposed three quantum integrable models for spin 1/2 chains with 
inhnite range interactions [l|. Two of these models, the rational or XXX and the hyperbolic- 
trigonometric or XXZ models, were solved for the spectrum and eigenstates. However, the 
exact solution of the third integrable model, the Elliptic Gaudin Model (EGM) or XYZ 
model, had to wait till 1996 2| for a complete solution in terms of the Algebraic Bethe 
Ansatz (ABA). The Gaudin models played an important role in the development and testing 
of several methods in quantum integrable theory, like the functional Bethe ansatz and the 

, ^ , the relation with the Knizhni^Zamolodchikov equations 


separation of variables 

n 

and with the corresponding Wess-Zumino-Witten models [7|, 
transformations js], etc. On a different perspective, the rational Gaudin model 


the construction of Backlund 

101 was 


linked to the exact solution of the reduced Bardeen-Gooper-Schriefer (BGS) Hamiltonian 
solved exactly by Richardson [u| and proved to be quantum integrable by Gambiaggio, 


Rivas and Saraceno 12| . Exploiting this link, several fami 
called Richardson-Gaudin (RG) models were proposed 


ies of exactly solvable models 


RG model found important applications to different areas of quantum many-boc 


including ultrasmall superconducting grains 


atomic gase s 


17 


ins QQ, 


10|. Since then, the rational 
y p hysics 
161, cold 


Tavis-Gummings models 


15 


19|, quantum dots [20|, l2l| and nuclear structure 


22 


(for a review 


see [24, l25[). More recently, the hyperbolic RG model was applied to describe p-wave 


pairing in 2D lattices 


26 


28| and ID Kitaev wires 29|]. With less success, there have been 


attempts to apply the EGM to matter-radiation problems including counter-rotating terms 


30, 


3l| . However, these integrable models lack of the radiative excitation term or lead to 


non-hermitian Hamiltonians. On a different respect, the EGM has been used to study the 


thermalization process o 


a spin chain with long range interactions in the transition from 


integrability to chaos 32|. None of these works attempted to hnd a numerical solution of 
the Bethe equations. 


The aim of this paper is to survey and generalize the exact solution of EGM for arbitrary 
SU(2) spins systems, and to study numerically the properties of the exact solutions in 
different scenarios. We hrst introduce the model with the corresponding exact solution in 
Sec. H. We then discuss a simple model of three different spins in Sec. HI for which we 
show how to solve the Bethe equations in order to obtain the complete set of eigenstates. 
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Next, we move to a physically oriented problem, a new integrable anisotropic central spin 
model (ACSM), in Sec. IV. In Sec. V we solve exactly the ACSM for long chains, and we 
extrapolate the ground state energy to the large N limit, showing that it coincides with the 
classical spin approximation in the thermodynamic limit. 


II. THE ELLIPTIC GAUDIN MODEL 

The EGM was hrst proposed by Gaudin|l| as a particular family of integrable spin Hamil¬ 
tonians with a fully anisotropic spin-spin interactions. The N commuting integrals of motion 
for a system of N arbitrary spins S'", with a = x,y,z and i = 1, ■ ■ ■ , N are 


N 


ID _ \ ^ Tx Qx Qx _i_ jy cy cy jz qz qz 


( 1 ) 


J = 1 

where the matrices satisfy the Bethe equations 

ja t7 I 7/3 77 I to t/3 _ q 

in order to fulhll the integrability conditions [i?j,i?j]=0. 

Following Ref. the can be expressed in terms of the doubly periodic elliptic Jacobi 
functions of modulus k, srL{z,k), crL{z,k) and dn{z,k) (for brevity, in general, we will not 
explicitly write the modulus), and a set of N arbitrary coefficients Zi as 


JZ = 




Jl = 


1 + k sn^{zi — Zj) 
sn{zi - Zj) 

1 — k sn^(zj — Zj) 


= 


*3 sn(zi - Zj) 

cn^Zi — Zj) dn{zi — Zj) 


( 2 ) 


sn(^i - Zj) 

Alternatively, the integrals of motion Ri can be expressed in terms of the raising and 
lowering spin operators dz iS^ as 


k 

- sn 


r, = Y. 

j(di) ^ 

cn^Zi — Zj) dn{zi — Zj) 


{z, - z,) (srs* + S-S-) + {srs- + S-S*) 


QZ QZ 

sn{zi-Zj) * \ 


( 3 ) 
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The elliptic integrals of motion break the su{2) symmetry associated with the conser¬ 
vation of the total spin S and the m( 1) symmetry associated with the conservation of the 
z component of the total spin S^. However, the model preserves a discrete Z 2 symmetry 
associated with a vr rotation of every spin around an arbitrary axis. Assuming 2 ; as the 
quantization axis, a rotation by an angle vr around this axis is related to the parity opera¬ 
tor P = Hill with eigenvalues -|-1 for positive parity and —1 for negative parity. 

Therefore, the complete set of common eigenstates of the integrals of motion can be classihed 
in two subsets of even or odd parity. 

The exact solution comprising the eigenvalues of the integrals of motion and the equations 
for determining the spectral parameters for a system of N spins with Sj = | has been 
obtained by means of the ABA in references 2, Here, we present the exact solution for a 
general system of N arbitrary spins s,. The derivation starting from a system of 1/2 spins 
is given in the Appendix. The eigenvalues of the integrals of motion ([T]), depending on a 
set of M roots to be determined by the Bethe equations that will be introduced below, 
are: 


A = Si 


N M 

Sj - Zj) + ipi{zi - Zj)) - {(p^{zi - K) + (pi{zi - A„)) i 


a=l 


ttI 


, (4) 


where M = combination of spins s* is allowed as long as the resulting summa¬ 

tion M is integer. Notice that M = N/2 for the Sj = 1/2 spin case and therefore, it would 
only admit an exact solution for an even number of spins N. In addition, I is an integer 
number that can take the values 0 or 1 in order to identify the two parity sectors. 


At this point, we have to introduce the Jacobi Theta functions [33|, [3^ 0 (m, k) = ■^4 {v\q) 
and H(m, k) = 'di{v\q), with the variable transformation v = ^. In these definitions K{k) = 
-j==== is the complete elliptic integral of the first kind, K'{k) = K (\/l — and 
the nome q = . The functions V5 i(m), and can be dehned now following 351, 


as : 


^i{u) = 


H'iu) 


ifi = 


0'l 


u] 


Q{u) 


The M roots Aq, in eq (j4]) are determined by solving the set of M coupled Bethe equations 
(see the Appendix): 
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N ^ 

Sj {^Pi{\a - Zj) + ^l{Xa - Zj)) - ((^4(Aa - Xfi) + (^l(A« - A/j)) + 1— = 0 . (5) 

i=i /3(^«) 

The function <y5(A) = <yCi(A) + </? 4 (A) appearing in the eigenvalues (jl]) and in the Bethe 
equations (j5]) has the special property of being periodic in the real part of its argument 
but “quasi periodic” in the imaginary part, i.e. Lp{X + i K') = </?(A) + zC(/c), where K' is 
the quasi-period (in the imaginary direction) and C{k) is a real constant only depending 
on the elliptic modulus k. As a consequence the imaginary part of the Bethe roots are 
constrained to its natural interval Im(A) G ^]- Numerical solutions obtained outside 

of this interval may lead to spurious non-physical results. On the other hand, as the real 
period of <p(A) is 2K, the Bethe roots of any physical solution should be constraint to the 
fundamental rectangle of the complex plain given by [0, 2K] x ■ 

We next analyze the behavior of the integrals Ri in the two limits A; —?■ 0 and /c —?■ 1. 
Taking into account that for /c —)■ 0, sn (u, k) —>■ sin (u), cn (u, k) cos (u) and dn (u, k) 1, 
it is easy to check that the elliptic integrals Ri ([3]) transform into the trigonometric ones: 




2 sinf^;,: 


— {S+S- + S-S+) + cot{zi - Zj) S!S] 
^j) 


On the other hand, starting from the x, y, z representation ([T]) and taking the limit /c —?■ 1, 
the elliptic functions transform according to sn (u, k) tanh (u) and cn (u, m ), dn(M, k) 
sech (u), thus we obtain: 


Di _ V_iI cosh^Zj - Zj) + sinh^Zj - Zj) ^ 

* ^ 2cosh(2;i — Zj) sinh(.Sj — Zj) ^ ^ ^ ' cosh(2;j — Zj) sinh(2;j — Zj) * ^ 

Performing a cyclic permutation of the axis and making use of some identities of the 
hyperbolic functions, we hnally obtain the Gaudin hyperbolic integrals: 


= E Uiuhfa-,,) + coth(,.- SfSJ , 

where rji = 2zi. 

In a similar way, it can be shown that the eigenvalues (j4]) and Bethe equations (jS]), reduce 
to the corresponding trigonometric and hyperbolic Gaudin solutions. 
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ZJX 

^12 

^12 

HI 2 

TJX 

^13 

uv 

^13 

Hh 

TJX 

^23 

uV 

-^23 


1.28522 

1.23563 

1.2293 

1.38861 

1.19509 

1.16865 

1.28522 

1.23563 

1.2293 


Table I. Exchange couplings of the three-spin Hamiltonian ([6|) . 


III. A THREE-SPIN SYSTEM 

In order to understand the behavior of the Bethe roots for the different eigenstates we 
treat in this section the simplest integrable EGM with three different spins, N = 3, and 
Si = 1/2, S 2 = 1, S 3 = 3/2. We construct an exactly solvable spin Hamiltonian as a linear 
combination of the integrals of motion 

3 3 

= ( 6 ) 

i<.j a=l 

We choose for the Hamiltonian H = —— jR 2 , with the parameters zt = 0, 0.2, 0.4 and 
the elliptic modulus k = \. The corresponding coefficients iJ" are displayed in Table [B 

The dimension of the Hilbert space is d = nil (2 Sj + 1 ) = 24. The Hamiltonian matrix 
is block diagonal with d_|_ = 12 for positive parity and d- = 12 for negative parity. The 
number of Bethe equations (|5]) and roots is M = X^li = 3. For each parity sector, 
dehned by / = 0 or 1 , we look for 12 different solutions with the three roots restricted to 
the fundamental rectangle defined by 0 < Re(A) < 2K and Im(A) < For this 

particular case 2K = 3.3715 and ^ = 1.07826. 

Tables [TTl and [TTIl depict the complete set of eigenvalues and the corresponding values of 
the Bethe roots for positive and negative parity respectively. For positive parity {I = 0) and 
real parameters Zi, the Bethe roots are real or complex conjugate pairs with the exception 
of roots having the imaginary part equal to half the imaginary quasi-period In the 

latter case the pair of complex roots need not to be a conjugate pair (real parts could be 
different). Moreover, complex conjugation of all roots gives the same solution due to the 
quasi-periodicity of the Jacobian functions. On the other hand, for negative parity (/ = 1) 
one of the roots is always complex with the imaginary part equal to half the imaginary 
quasi-period (+^), to compensate the imaginary term added to the Bethe equations. The 
other roots could be real, or complex pairs following the same behavior as in the I = 0 case. 

As a hnal remark, we have checked that the eigenvalues obtained solving the Bethe 
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E 


^2 

^3 

-8.13147 

0.277673 

0.261164 -h 0.115827 i 

0.261164 - 0.115827 i 

-5.64950 

0.0674548 

0.268200 

2.150100 

-5.48168 

1.92346 

0.281146 -h 0.0491301 i 

0.281146 - 0.0491301 i 

-0.850805 

0.286691 

0.726263 1.07826 i 

3.15855 - 1.07826 i 

-0.758290 

0.286760 

0.257994 1.07826 i 

1.94100 - 1.07826 i 

-0.649792 

0.0448926 

2.06330 -h 0.519302 i 

2.06330 ~ 0.519302 i 

-0.615993 

0.0473203 

0.375176 1.07826 i 

2.06325 - 1.07826 i 

-0.606659 

0.0481756 

0.825910 - 1.07826 i 

3.29741 -h 1.07826 i 

-0.394121 

0.287018 

1.94224 -h 0.516064 i 

1.94224 - 0.516064 i 

6.69616 

1.95117 

0.612035 - 1.07826 i 

3.29405 -h 1.07826 i 

6.88050 

1.95157 

0.268084 - 1.07826 i 

1.95185 -h 1.07826 i 

7.22436 

1.95228 

1.95248 -h 0.704049 i 

1.95248 - 0.704049 i 


Table II. Eigenvalues and Bethe roots in the positive parity sector {I = 0) of the three-spin Hamil¬ 
tonian dH). 

equations ([5]) fully agree with the results of an exact diagonalization of the Hamiltonian (j6]). 


IV. THE ANISOTROPIC CENTRAL SPIN MODEL 


The central spin model (CSM), describing the hyperhne interaction of an electron spin 
in a quantum dot with a non-interacting system of surrounding nuclear spins, has been 


proposed as the main component of spintronic devises and solid state qubits 


3^. The 


isotropic CSM Hamiltonian with Heisenberg exchange couplings between the central spin 
and the nuclear spin bath, subjected to an external magnetic held, is precisely one of the 
integrals of motion of the rational RG model. As such, it has been extensively studied using 


exact solutions 




Anisotropic effects due to the hyperhne interaction between the 


central spin and the nuclear bath can still be investigated within the Hyperbolic or XXZ 
RG model j^. However, the inclusion of the quadrupole coupling in the electron-bath 


interaction goes beyond the XXZ model 


Here, as a physically oriented example of a 
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E 


^2 

^3 

-5.86850 

0.239016 -h 1.07826 i 

0.280492 -h 0.0487235 i 

0.280492 - 0.0487235 i 

-5.64331 

0.0682364 

0.268316 

2.14920 -h 1.07826 i 

-5.59109 

0.0727121 

0.269232 

0.458056 -h 1.07826 i 

-5.52799 

0.281070 -h 0.0490841 i 

0.281070 ~ 0.0490841 i 

1.92361 1.07826 i 

-0.714459 

0.286793 

1.94080 

0.258161 -h 1.07826 i 

-0.649689 

0.0448967 

2.06348 

2.06312 1.07826 i 

-0.619842 

0.0469765 

2.06361 

0.375167 -h 1.07826 i 

-0.395533 

0.287017 

1.94206 

1.94243 1.07826 i 

6.59546 

0.267512 -h 1.07826 i 

0.985596 -h 1.07826 

2.91839 - 1.07826 i 

6.64346 

0.599786 1.07826 i 

1.95133 - 1.07826 i 

3.30613 1.07826 i 

6.88448 

0.268092 -h 1.07826 i 

1.95170 -h 0.496275 i 

1.95170 - 0.496275 i 

7.22431 

1.95235 -h 0.345399 i 

1.95235 ~ 0.345399 i 

1.95255 -h.07826 i 


Table III. Eigenvalues and Bethe roots in the negative parity sector (I = 1) of the three-spin 
Hamiltonian ([61). 


quantum integrable system derived from the EGM, we study a modified anisotropic central 
spin model (ACSM) without an external magnetic field. The introduction of an external 
magnetic held breaks the integrability of the EGM since it does not admit linear terms as 
opposed to the rational and trigonometric-hyperbolic cases. In our model the system is 
composed by a single electronic spin si = | interacting with the N — 1 nuclear spins | of the 
bath. The hyperhne and quadrupole couplings between the electron and the surrounding 
spins is described by a completely anisotropic antiferromagnetic exchange interaction. We 
assume that the electron spin is located at position zi = 0, while the environmental spins 
are uniformly distributed within a linear segment at a hnite distance a {z 2 = a) with the 
last spin at position z^ = b, with 0 < a < b < K. The restriction for b to be lower than 
or equal to K is necessary to keep the interaction decreasing with distance in the selected 
interval. Therefore, the values of the fixed parameters 2 ; are given hy zt = a + j^{b — a) 
for i > 1. The anisotropic central spin Hamiltonian is defined by the first integral of motion 
(d]) of the EGM: 




















Figure 1. Exchange couplings of central spin Hamiltonian ([7]) as a function of distance. In the 
horizontal axis we show as an example a system with N = 8 spins. 


N 

Hacsm = -Ri = Y^ + Jp\s] + r^SlS^) , ( 7 ) 

J=2 

where J" = = Jij{zj) for a = x,y, z as given in Eq. ([2]). The properties of elliptic 

functions determine the relation between the exchange interactions in the x, y, z directions 
as Jf > > J? for 0 < Zj < K. Fignre [U shows the three components of the interaction 

as a fnnction of distance for k = 0.5. In the horizontal axis we display, as an example, a 
spin network with iV — 1 = 7 environmental spins nniformly distribnted in a segment with 
a = 0.2 and b = 0.6. 

In order to gain insight into the strnctnre of the GS wavefnnction we explore the clas¬ 
sical description of the model. In this approximation each spin ^ operator is replaced by 
S = ^[cos(p sin9, sinip sin9, cos9], where 9 and (p are the nsnal polar and azimnthal angles 
satisfying 0 < 6* < vr and 0 < < 27r. Inserting into the Hamiltonian ([7]), the classical 
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energy is given by 



( 8 ) 


Our task is to find the absolute minimum of the classical energy E^i with respect to the 
angular variables of the central spin 6*i, ipi, and those of the nuclear spins {9j, Each 


term of the sum Ej depends exclusively on the angles of a particular nuclear spin j and 
the angles of the central spin. Minimization with respect to the angles leads to a system of 


2 X N coupled variational equations. However, minimizing each term Ej independently will 
yield the absolute minimum if each solution is compatible with a unique value of the central 
spin variables Oi,^pi. Notice that even though the spins in the bath are non-interacting, 
the latter condition induces an effective interaction through the central spin. By solving 
the four variational equations derived from the minimization of each Ej we obtain different 
classes of solutions corresponding to stationary values of Ej in the set ^Jj^ }• 

As the largest component of the interaction is Jj we conclude that the minimum for each 
term is Ej = —Jj, and the corresponding angles are 9i = |, (pi = vr, 9j = |, (pj = 0 (j > 1). 
The classical GS corresponds to an antiferromagnetic state with all spins aligned into the 
X axis, and the central spin pointing in the opposite direction to the environmental spins. 
Therefore, for the minimal energy configuration, the classical energy per spin is: 




In order to find an expression for the classical energy density in the thermodynamic 
limit, i.e. limAr^oo "yf, we define a normalized density distribution for the spins such that 


p(z) dz = 1, where hi C (0, K) is the compact interval containing all environmental spins 
(the parameters z^s except for zi), so that the number of spins in an interval of length dz 


is given by (iV — 1) p{z) dz. Introducing this distribution in equation (j9]) and taking the 
corresponding limit we obtain in general: 



( 10 ) 


According to our model of equidistant bath spins, we assume a uniform distribution in 
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the interval hi = [a, b]: 


p(^) = 1-^ ■ (11) 

b — a 


( 11 ) 


Finally, from Eq. ffTOj) and making use of the uniform density flTT]) we obtain for the classical 
energy density in the thermodynamic limit: 



( 12 ) 


Q 


which can be integrated to give 


1 , r sn(b) (cn(a) + dn(a)) (dn(&) — k cn(&)) 

- l02( - 

4(a — b) [sn(a) (cn(6) + dn(6)) (dn(a) — k cn(a)) 


(13) 


We will later compare this classical energy density with the large N limit of the exact 
solution. 

V. EXACT SOLUTION OF THE ACSM 

Let us now turn our attention to the exact quantum solution of the ACSM. We will 
analyze the distribution of spectral parameters or Bethe roots Aq of the Bethe equations 
(jS]), as well as the ground state energy (per spin) of the Hamiltoninan Hacsm, for several 
system sizes up to A = 300. A large N extrapolation of these results will allow us to compare 
with the classical energy density. We choose the modulus k = 0.5, thus the fundamental 
interval is defined by iF(0.5) ~ 1.68575. For the 2 ;’s interval we set a = 0.2 and b = 0.6. 

The numerical solution of the nonlinear system of equations ([5]) faces the usual problems 
of any Gaudin system, namely a huge number of independent solutions (the dimension of 
the system grows exponentially with N) and dangerous divergences which hinders numer¬ 
ical procedures based on iterative methods. Moreover, finding a specihc solution strongly 
depends on the initial guess. To overcome these difficulties we first solve the system for a 
small number of spins, where we can identify the distribution of roots for the ground state 
as well as for all the excited states. We then follow a particular solution increasing N by 
means of an algorithm described below. 

In Fig. [2] we show the distribution of Bethe roots in the complex plane for the ground 
states (/ = 0) of two systems, a small system with A = 12 and large system with A = 300. 
The complex plane is divided into three regions delimited by vertical lines crossing the real 
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Figure 2. Ground state Bethe roots of the ACSM ([Tj) for different system sizes. 

axes at Re{\) = 0, a, b and 2K. We found that the ground state has the root Ai, which 
is associated with the central spin, always real and located in the hrst region satisfying 
0 < Ai < a, while the other roots are located over a smooth curve, symmetric with respect 
to the real axes, with b < Re{X) < 2K. We can see in Fig. [2l for N = 12, that the hve 
roots are almost vertically aligned at Re{X) ~ 2.120. For N = 300 the 149 roots make up 
an almost vertical segment at Re{X) ~ 2.087. The inset amplihes the third region adding 
several solutions for the intermediate systems with with N = 20, 40 and 100. For increasing 
values of N the real part of the roots in the arc decreases approaching a limiting value. 
Moreover, the hrst root Ai —)■ 0 in this limit as can be seen in Table IIVI 

The ground state always displays the same pattern, with one real root in the hrst region 
0 < Ai < a and iV — 1 roots distributed over a smooth arc in the third region, i.e. Re{X) > b. 
A similar pattern takes place for the lowest energy state in the I = 1 block. In Fig. [3] we 
show the distribution of Bethe roots for the GS (already displayed in Fig. [2]) and the hrst 
three excited states in positive parity sector of the ACSM with 12 spins. As it can be seen, 
the arc of complex roots which characterizes the GS persists for the low lying excited states 
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N 

XiN 

Min(Re(A)) 

Max(Re(A)) 

E 

N 

12 

0.327300 

2.120138 

2.120487 

-0.821616 

20 

0.329466 

2.106036 

2.106256 

-0.792253 

40 

0.330390 

2.095787 

2.095898 

-0.772954 

80 

0.330681 

2.090745 

2.090800 

-0.7640486 

100 

0.330726 

2.089743 

2.089787 

-0.762322 

200 

0.330805 

2.087743 

2.087765 

-0.758920 

300 

0.330828 

2.0870780 

2.0870926 

-0.757801 

oo 

0.330869 

2.0857505 

2.0857505 

-0.755586 


Table IV. Extremum real values of the GS state Bethe roots and corresponding energies per spin 
for different N values. 


'b 

HH 



0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 


Re(A,) 


Figure 3. Bethe roots of the GS and hrst 3 excited states of the ACSM with N=12 spins. 

while some detached roots are distributed in other regions of the complex plane. Different 
distributions of the roots give rise to the complete set of eigenstates. 

The scale of the graph in figure El does not allow to appreciate the detailed form of the 
arcs in the third region. In hgure |4] we show, with a smaller scale, the actual shape of the 
arc for N = 300. The difference between the maximum and minimum real parts of the 
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Figure 4. Ground state Bethe root of the ACSM d?]) for N = 300. The continuous line corresponds 
to the fit of Eq. (fTT)l . 

roots in the arc decreases with increasing N, becoming null at —?■ oo. Therefore, in the 
continuous limit the arc turns into a vertical segment with Re(A) = 2.08575 with half of the 
quasi-periods as the interval extremes < Im (A) < (see table HV|) . 

In order to obtain the numerical ground state solution for a large N system, we start 
from the solution of a small system. We then increase N, typically doubling it in each 
iteration. In each step we make a least square ht of the complex arc expressing the real part 
Re(A) = X as a function of the imaginary part Im(A) = y taken as the independent variable. 
An excellent agreement is obtained for any N value by means of the 4-parameter function: 

X = a + (3 dn(ci y) cn(c 2 y) . (14) 

The continuous line in Fig. 0] corresponds to a ht of this function for the N = 300 system. 
We take advantage of this excellent ht to generate, for each system of size N the initial guess 
from the lower size system Nq. In addition, the initial guess for the hrst root Ai is obtained 
by a linear scaling Ai = A?^. In both cases the index 0 stands for the solution of the 
previous step. The procedure is stable, and it allows to hnd numerical solutions for very 
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N 


Figure 5. Exact quantum and classical GS energies per spin of the ACSM for different N values. 
The extrapolated —)• oo limits coincide with the analytic value (|13ll . 

large systems. 

It is interesting to analyze the numerical results in the N ^ oo limit. Table HV] shows 
these results for several values of N and the numerical extrapolation to the thermodynamic 
limit. The hfth column of the table displays the energy per site and the extrapolated value 
in the thermodynamic limit. We can now compare these results with the classical energy 
density. For a = 0.2, b = 0.6 and k = 0.5 equation ffTOll yields Sd = —0.75558603 showing 
an excellent agreement with the extrapolated result. In hgure ([5]) we show a comparison 
between the classical energy ([9]) (open circles) and the quantum energy (black circles) for 
several values of N. Both branches converge to a common limit for iV —)■ oo. Continuous lines 
correspond to a third order polynomial least square hts. The gap between the classical and 
the quantum energies for hnite systems is a direct consequence of the quantum fluctuations 
that disappear in the thermodynamic limit. 


VI. SUMMARY 

The rational and hyperbolic Gaudin models have been extensively employed lately to 
study many-body quantum systems in several branches of mesoscopic physics. On the 


15 



contrary, the EGM derived by Gaudin in 1976 and solved exactly with the ABA in 1996 
has been scarcely used as a mathematical tool to investigate many-body physical problems. 
In this article we have hrst summarized the integrals of motion of the model and the exact 
solution for systems with arbitrary spins. Subsequently, we discussed the behavior of the 
Bethe roots for the complete set of eigenstates of a system with three different spins. We 
showed that the Bethe roots should be restricted to the fundamental rectangle in order to 
warrant that every solution corresponds to a physical eigenstate. Finally, we introduced a 
particular anisotropic GSM, accounting for the hyperhne interaction of an electronic spin in 
a quantum dot with the environmental nuclear spins, and possible effects due to quadrupole 
couplings. The so called AGSM was solved exactly for large number of spins and the GS 
energy was extrapolated to the thermodynamic limit, which coincides with the classical 
approximation. We hope that our numerical studies would pave the way to applications of 
the EGM to other physical many-body systems. 
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APPENDIX. GENERALIZATION OF THE EGM TO ARBITRARY SPINS 


We start with the known inte, 
equations for spin 1/2 systems 


^als of motion, the corresponding eigenvalues and the Bethe 

6 | 


N 


{Za - z,) sisi + r {Za - z,) sisi ], 

b(^a)=l 


ha = ^ ^ {Za - Zb) - - {Za - K) + 

h(^a) = \ a.=l 


4a:’ 


^ M M . j 

Za) - + ^ = 0’ 

a=l /3(^a)=l 

where M is the number of spins and M. = Mj‘2 the number of roots. 


(16) 


(16) 


( 17 ) 
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Following Ref. [l| we now group an arbitrary number of adjacent spins into clusters and 
define a new lattice with sites i and spins Si = YTa=ai where a* is leftmost site of the 
cluster i and n* the number of spins in the cluster. The new lattice fulhlls Af = 
where N is the number of clusters or sites in the new lattice. We recall that the functions 
J {x) and (fi (x) are odd functions. The new integrals of motion are 


ai+rii ai+rii ai+Ui 

R, = Ha = {Za - z,) SlSl + r {Za - Z,) S^Sl] + 

CL—di CL—CLi 

CLi+rii 

Y, E1-^“' + J" • 

a=ai b^Ci 

where Ci = {oj, Oj+i, • • • , Oj + n*}. The first term in the right hand side cancels out due to 
the antisymmetry of the functions J (x). We now assume that the parameters z inside each 
cluster are all equal, Za = Zi for all a & Ci, and together with the definition of the cluster 
spins Si we obtain the integrals of motion in the general case 

N 

H - ^3) - ^3) + J' - ^3) S!S]] . ( 18 ) 

3 {¥=i )=3 

The corresponding eigenvalues and Bethe equations are transformed as 

ai+Ui N M , 

c = ^ ha = ^ rij^p (Zi -Zj)-^Y^ (19) 

a=ai j(:^i)=l a=l 

and 


1 ttI 

-Y^i^i^Oi-Zi)- Y (^« - ^/ 3 ) + = 0 . ( 20 ) 

i=i 0{^C=1 

It can be shown pj that the exact solution corresponds to the maximum spin in each 
cluster, namely s* = y, which implies that M = Ai With these last 

replacements we obtain the eigenvalues (jlj) and Bethe equations fl5l) . 
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